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ABSTRACT 

We present an algorithm for selecting gravitational lens candidates from 
amongst Sloan Digital Sky Survey (SDSS) quasars. In median Early Data Re- 
lease (EDR) conditions, the algorithm allows for the recovery of pairs of equal 
flux point sources down to separations of ~ 0'.'7 or with flux ratios up to ~ 10:1 at 
a separation of 1'.'5. The algorithm also recovers a wide variety of plausible quad 
geometries. We also present a method for determining the selection function of 
this algorithm through the use of simulated SDSS images and introduce a method 
for calibrating our simulated images through truth-testing with real SDSS data. 
Finally, we apply our algorithm and selection function to SDSS quasars from 
the EDR to get an upper bound on the lensing fraction. We find 13 candidates 
among 5120 z > 0.6 SDSS quasars, implying an observed lensing fraction of not 
more than 4 xlO~^. There is one likely lens system in our final sample, implying 
an observed lensing fraction of not less than 3 x 10~^ (95% confidence levels). 

Subject headings: quasars: general; gravitational lensing 

1. Introduction 

The discovery of new gravitationally-lensed quasars, where the gravitational potential 
of a foreground galaxy produces multiple images of a background quasar, is of broad as- 
trophysical interest. Such systems can be used individually to directly model the mass 
distribution of the lens galaxy (Kochanek 1995; Cohn et al. 2001; Munoz, Kochanek, & 
Keeton 2001; Koopmans & Treu 2002; Rusin et al. 2002; Treu & Koopmans 2002; Trott & 
Webster 2002 ), investigate dark matter subtructure (Metcalf & Madau 2001; Chiba 2002; 
Dalai & Kochanek 2002), measure the Hubble constant through time delay measurements 
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(Refsdal 1964; Kundic et al. 1997; Koopmans & Fassnacht 1999) and to probe the source 
size through microlensing studies (Wyithe et al. 2000). Additionally, the statistics of an ob- 
served sample of gravitational lenses can be used to infer the mass profiles of the population 
of lens galaxies (Keeton 2001; Rusin & Ma 2001; Takahashi & Chiba 2001) or to constrain 
the value of cosmological parameters, in particular the value of the cosmological constant 
(Turner 1990; Fukugita, Futamase, & Kasai 1990; Kochanek 1996; Helbig 1999; Chae et al. 
2002; although see Keeton 2002). This last apphcation is of particular relevance in hght 
of the recent evidence of a high value of A from measurement of high-redshift supernovae 
(Perlmutter et al. 1999; Riess et al. 2001), CMB anisotropy measurements (Netterfield et al. 
2002; Halverson et al. 2002; Lee et al. 2001) and so-called concordance models (Bahcall et al. 
1999). Lensing statistics will almost certainly never be competitive with the accuracy of the 
aforementioned experiments in terms of measuring the values of cosmological parameters. 
Nevertheless, confirmation of the high lensing rates predicted in high A model cosmologies is 
an important observational test, particularly as the systematic errors for lensing are distinct 
from those of supernovae experiments. 

The Sloan Digital Sky Survey (York et al. 2000, henceforth SDSS) datasct offers an 
opportunity to discover a large number of new gravitationally-lensed quasars. Primarily, it 
is the size of the survey which makes it a superlative dataset for finding lenses. Assuming a 
plausible lensing fraction of ~ 10~^, the complete spectroscopic sample of 10^ SDSS quasars 
will include of order 100 lens systems. Further, the 10^ square degrees of imaging data should 
contain of order 10^ quasars within one magnitude of the spectroscopic limit, although their 
lensing rate will be probably be considerably lower, depending on the details of the quasar 
luminosity function. A number of other features of the survey are also pertinent. Since 
the imaging data arc taken in five bands nearly simultaneously, it is possible to use colour 
information and require consistency between different bands to considerably improve the 
efficiency of most selection techniques for lenses. Further, the SDSS is particularly well- 
suited to the collection of a statistically well-defined sample of lenses because the uniformity 
of the data and of the reductions should allow for an accurate determination of the selection 
function. In this paper we present the methods by which we intend to assemble such a 
sample as well as the techniques we will use to measure our selection function through the 
use of simulated images. At present, we consider only the problem of recovering lenses from 
a sample of objects which are already known to contain at least one quasar component from 
their SDSS spectra. We will not explicitly treat the problem of how the SDSS quasar target 
selection may influence the selection function (Kochanek 1991, henceforth K91), although 
we do discuss this problem in the following section. 

When discussing lensing statistics, we deflne lens fraction to mean the true, unknown, 
fraction of quasars in the sample which are multiply-imaged and we define observed lens 
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fraction to mean the fraction of quasars in the sample which have been observed to be 
multiply-imaged. The observed lensing fraction is related to the lensing fraction through the 
selection function. Note that, even with complete knowledge of the selection function, the 
lensing fraction cannot be recovered from the observed lensing fraction as there may always 
be a glut of unobserved lenses in the selected-against regions of parameter space. However, 
given a theoretical model of the lensing population and cosmology, one can use the selection 
function to predict the observed fraction and compare to the experimental value. 

Readers unfamiliar with the SDSS are encouraged to consult the SDSS EDR paper 
(Stoughton et al. 2002) in case of encountering unfamiliar SDSS terminology. 

The structure of this paper is as follows. In section 2 we briefly consider the problem 
of selecting Icnscd quasars in general and discuss those aspects particular to the SDSS. We 
then introduce a selection algorithm based upon fitting the observed point spread function 
(henceforth PSF) to the SDSS atlas images, and illustrate how the selection function of this 
algorithm can be measured. In section 3 we apply our selection technique to a sample of 
SDSS quasars to show how a sample of lens candidates for high-resolution follow-up can be 
defined. In section 4 we discuss our results and consider future work. 



2. Selection Techniques and Selection Functions 

2.1. General Considerations 

A sample of SDSS lenses will be of greatest interest if it can be selected by a method 
which could be applied uniformly to the entire SDSS dataset. In an ideal lens survey, every 
SDSS quasar would be foUowed-up with high-resolution imaging to test for the presence of 
multiple images. However, given that the complete spectroscopic sample of SDSS quasars 
is scheduled to comprise some 10^ objects, such an approach is not immediately feasible. 
An obvious alternative is to use available morphological information from the SDSS imaging 
data when selecting candidates for follow-up (Webster et al. 1988). In this way the com- 
pleteness of the survey is reduced, but the efficiency is expected to be greatly increased. 
Before describing the mechanics of our morphological criterion, we briefly review some of 
the theoretical expectations and observational constraints which motivate our approach. 

The precise details of the predicted distribution of lens image splittings depend upon 
the mass model and luminosity function chosen for the population of lens galaxies. However, 
it is a generic prediction of reahstic models that the peak of the distribution is at about 1" 
and that the majority of systems will exhibit sphttings less than 2" (Turner, Ostriker, & 
Gott 1984, Hinshaw & Krauss 1987). Results from the CLASS radio lens survey are in good 
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general agreement with such expectations (Browne et al. 2002). Yet the typical FWHM for 
SDSS Early Date Release (Stoughton et al. 2002) imaging data is not better than I'.'G and 
since the SDSS pixel scale is 0'.'4, critical sampling requires that seeing not be significantly 
sub-arcsecond. Consequently, to recover a large fraction of lens systems, it is necessary to 
be able to identify unresolved systems. Although the SDSS photometric pipeline (Lupton 
et al. 2001, henceforth PHOTO) does include software which deblends compound objects, 
we find that few pairs with splittings below 3" are deblended automatically by the pipeline. 
In this work we concentrate on methods for recovering such small-splitting lenses. 

For wider separation lenses (A^^ > 3"), we expect components to begin appearing as 
separate objects in the SDSS catalogs, and a search for such lenses can be performed by 
matching pairs in the quasar target lists. It should be noted that fibre collisions prevent the 
SDSS from obtaining separate spectra on a single plate for pairs of objects with a separation 
of less than 55". There is a considerable amount of overlapping area between different plates 
(~40%), so that two nearby objects could be observed on separate plates, but many such 
close targeted quasar pairs will have only one SDSS spectrum, meaning that even at these 
separations most lenses cannot be identified from SDSS data alone. 

In addition to separation, the multiplicity and configuration of lens images also greatly 
affects their detectibility. Almost half (10 of 22) of the lenses found by the CLASS survey 
exhibit four or more images. Ideally, our selection techniques would be effective at finding 
lenses of any configuration. However, there are two main reasons why this is presently not the 
case; i) The computational load associated with techniques which might be best suited for 
finding higher multiplicity systems (ie additional model components, deconvolutions) begins 
to become prohibitive for the thousands of objects we wish to analyze. Our morphological 
criterion (described below) is undoubtedly more sensitive to two than four image systems, ii) 
The modeling of quads is obviously more involved and less certain than the modeling of two- 
image systems. In fact, Rusin & Tegmark (2001) argue that simple elliptical models cannot 
reproduce the quad fraction observed by CLASS. As a result, the selection function for quads 
is difficult to characterize. In section 2, we explicitly calculate the selection function for pairs 
of point-sources but we are only able to roughly estimate the detectibility of quads. 

In the observational section of this paper, we will only consider selecting lenses from a 
parent sample of objects which have already been confirmed to contain at least one quasar 
component by SDSS spectroscopy. This is largely because in our early efforts to select lenses 
from the photometeric data alone we found a high number of false positives, in particular 
pairs of blue compact emission-line galaxies which are quite difficult to distinguish from 
quasars on the basis of broad-band colours alone. However, having effectively restricted 
ourselves to the SDSS quasar sample, an important question, as pointed out in K91, is 
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whether the selection of the parent quasar sample is biased with respect to lenses. Even 
if we only consider the quasar images themselves, we already find that the answer is yes. 
The SDSS quasar target selection (Richards et al. 2002, henceforth R02) effectively requires ^ 
that quasars with a redshift greater than 3.0 be classified as point sources by PHOTO. While 
this does not eliminate all recoverable lenses, as will be shown in the next section, it does 
represent a substantial bias. Additional biases can be introduced by the lensing galaxy. The 
lensing galaxy can cause a lens system to be excluded from the quasar sample either because 
i) the observed colours of the lens system are changed by flux from the lensing galaxy, ii) 
the observed shape of the lens system is changed by flux from the lensing galaxy, or iii) the 
observed colours of the lens system are changed by reddening of the quasar images along 
different sight-lines in the lensing galaxy. 

For most lenses, we expect the lensing galaxy to lie at moderate redshift. Turner et al. 
(1984) find that the peak of the redshift distribution of lenses hes at z ~ 0.5, with the 
majority of lenses lying at redshifts greater than 0.5. For high A cosmologies, the lenses are 
expected to lie at moderately higher redshifts (cf. Fukugita et al. 1990) This result is fairly 
model-independent since the term involving the angular-diameter distances in the expression 
for the critical surface density 



^cr OC ^ ^ (1) 

(e.g. Scheider, Ehlers, & Faico 1992) always tends to make lensing most likely when 
the lens lies roughly mid-way between the source and the observer. We can get a feel for 
the contribution of such a lens if we consider that the expected SDSS i magnitude of an 
galaxy z — 0.5 is 19.6 (20.6 and 22.3 in r and g, respectively), using Mg* — —20.04 
(Blanton et al. 2001), X-corrections and colours from Eisenstein et al. (2001), and including 
luminosity evolution in the amount — —0.6 (e.g. van Dokkum et al. 2001) with 

Q,M = 0.3, Qa = 0.7 and Hq = 75 km s~^ Mpc~^ For low redshift quasars {z < 3) SDSS 
quasar target selection requires that the i magnitude of an object be less than 19.1 in order 
to be included in the spectroscopic sample (R02), hence we expect the lensing galaxy to 
contribute significant flux to the SDSS atlas image in i, but less so in the bluer bands. In 
particular, the u — g colour, the chief identifier of low-redshift quasars, of the low-redshift 
sample will not be significantly affected by lens galaxy flux. However, for the fainter, higher- 
redshift sample [z > ?)) the i magnitude hmit is 20.2, and the objects are chiefly detected in 
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i and z, so the flux from the lens galaxy will likely make a substantial contribution. On the 
other hand, lensed quasars at higher redshifts will tend to have higher redshift, hence less 
prominent, lensing galaxies. It should also be noted that the above galaxy magnitudes are in 
Petrosian magnitudes, while quasar target selection uses PSF magnitudes. This means that 
target selection will preferentially consider concentrated quasar light over extended galaxy 
light when making colour cuts, reducing contamination from the lens galaxy. Also, quasar 
target selection does not simply select objects with 'normal' quasar colours, but rather targets 
objects which are outliers from the stellar locus. Hence, in order to cause a lens system to 
be de-selected, flux from the lens galaxy must not only change the colours, but must also 
change them in such a way that the object will move onto the stellar locus (or into one of the 
other excluded areas, see R02). The extent to which this is possible depends largely upon 
the redshifts of the quasar and the lens galaxy (cf. Figure 4 of R02 to see the location of the 
galaxy locus in SDSS colour space). 

A similar situation exists with regard to the effect of the lens galaxy upon the observed 
morphology. For the brighter, low-redshift sample, lens galaxy flux will not remove an 
object from the sample as no morphological cut is made in quasar target selection. For 
higher redshift objects, only unresolved objects will enter the sample, so systems where the 
lensing galaxy is prominent will again be excluded. 

Corrections for both of these selection effects must be made before any sample of lenses 
chosen from the quasar sample can be used for statistical purposes. We would argue that 
the use of simulated images as described in the next section is the best way to determine 
the required corrections. However, in this work we will not address this problem explicitly 
and henceforth we make the simplifying assumption that the flux from the lens galaxy 
is negligible. Of course, there is a particularily interesting class of lenses for which this 
assumption will fail quite badly, perhaps most notably for 2237+0305 (Huchra et al. 1985). 
However, the problem of finding lens systems for which the lens galaxy lies at low redshift 
is one which is possibly best approached through the use of the SDSS main galaxy sample 
(Mortlock & Webster 2000). 

Falco et al. (1999, henceforth F99) studied the effects of reddening in known lens systems 
and estimated the mean extinction to be E[B — V) = 0.10 ± 0.08 mag. This amount of 
reddening is comparable to the width of the low-redshift quasar locus in color space (cf. 
figure 5 of Richards et al. 2001), and so should not remove a large fraction of lenses from 
the quasar sample. In section 3, we present a selection technique which requires consistent 
colours between lens candidate components while allowing for differential reddening. 

In summary, none of the effects associated with the lens galaxy are predicted to remove 
a large fraction of lenses at redshifts less than 3, where quasars will be significantly brighter 
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than the lens galaxies, at least in the bluer bands, and where quasar target selection is 
unbiased with respect to object morphology. For the higher-redshift sample, flux from the 
lens galaxy will produce a considerable anti-lens bias both due to colour contamination and 
altered morphology. 



In this and following sections we present results based largely upon simulated SDSS 
images. The details of how these simulations are produced are found in the appendix. 

We devised a selection technique for selecting non-point source morphologies based upon 
fitting the observed PSF to the SDSS atlas images. Specifically, we fit each object twice; 
once with a model consisting a single PSF, and then again using a model consisting of two 
independent PSFs. Each PSF is allowed to vary in position and flux (normalization), so 

that the first model has three parameters and the second has six. The fits are performed 
using the empirical PSF as measured by the SDSS photometric pipeline. More accurately, 
since the photometric pipeline measures PSF variation across the frame, the PSF used is one 
re-constructed according to the object's position (see Lupton et al. 2001 for details). For 
each of these models, we minimize the reduced chi-squared: 



where Xi arc the object counts, Ui are the model counts, sky is the measured sky value, 
dark is measured dark variance, and iV is a sum over pixels in the atlas image. Note 
that the counts in SDSS atlas images are less the measured sky value. The minimizations 
are performed using the downhill simplex method (Press et al. 1992). We seed the two- 
component fit by measuring the inverse-radius weighted moments (1/r^) (cf. R. H. Lupton, 
in preparation) of the atlas image and taking as our initial guess an equal flux pair with the 
same moments. We can also measure the signal-to-noise (S/N) for an object as 



2.2. Morphological Selection 




N 



{xi + sky) I gain + dark 



(2) 



N 




(3) 



By summing over the entire atlas image we do not get optimal signal-to-noise (in the 
sense of, for instance, Naylor 1998), however, our measure has the advantage of being inde- 
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pendent of whether we are fitting a one or two component model. For simulated data the 
fits presented are in a single band, whereas for the real SDSS atlas images we performed fits 
independently in r and i, and, unless otherwise stated, the results of the fits are presented 
as averages over these three bands weighted by the square of the S/N. 

Our selection algorithm is based upon the difference in the residuals of our one and two 
component fits. Specifically, if we define and x| to bethe best-fit value of for the 1 or 2 
component model, respectively, then our morphological criterion depends upon the statistic 
Ax^ where 

AX^ ^Xl-Xl (4) 

Clearly, for a noise-free point source, we expect Ax^ = 0. We chose a threshold value 
of Ax^ which we determined would avoid any false positives from actual point source due to 
Poisson and sky noise. This value was obtained by fitting to, and hence, measuring the fluxes 
of, 5120 z>0.6 SDSS quasars and then simulating a sample of 5120 equivalent point sources 
with the same fluxes, seeing, and observed sky values. Figure 1 shows the distribution of 
Ax^ for this simulated sample. A Gaussian flt to this distribution predicts that even for the 
complete sample of 10^ quasars, the expected number of false positives is less than one for 
Ax^ > 0.05, although it is clear that such a fit does somewhat underestimate the number of 
objects at larger values of Ax^- We conservatively set our threshold value to be Ax^ > 0.1, 
and expect no false positives due to photon and sky noise. It should also be noted that of 
our simulated 5120 objects, 4 did show Ax^ > 0.1, but these events were associated with 
numerical issues in the fit and converged to lower Ax^ if the simplex fits were restarted. 

2.3. Characterizing the Selection Function 

2.3.1. Two- Component Systems 

In order to characterize the selection function which results from the selection criterion 
described above, we simulated an ensemble of SDSS images with two point-source objects of 
varying separation and flux ratio. We sampled the separation direction every 0.2 arcseconds 
between separations of 0.2 and 3.2 arcseconds and we sampled the flux ratio direction every 
factor of 2 between 1 and 64. At every lattice point in the separation - flux ratio parameter 
space we produced ten realizations of the geometry with random orientations, displacements 
relative to the pixel grid, and realizations of the noise. We also folded-in variations in the 
observed sky in the ten iterations (see appendix). For our fiducial set of simulations, we 
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assigned a total flux equal to the median total flux in r observed in our sample of 5120 
SDSS quasars (corresponding to a psf magnitude of r ~ 19.0) and we set the seeing to 
be equal to the observed median seeing (1'.'6 arcseconds FWHM, also in r) for this same 
sample. Note that when changing the flux ratio of a pair of point sources we conserved the 
total flux between them such that all simulated objects would appear identical if completely 
unresolved. 

Figure 2 shows the regions of the separation - flux ratio parameter space in which objects 
pass our selection criterion. Note that failure to select a object can cither result from the 
Ax^ in fact being less than 0.1 for the correct geometry, or from the two-component fit failing 
to find the correct geometry and corresponding global maximum of Ax^ ■ This is the case at 
larger separations and flux ratios (eg 3'.'0 and 4:1). Such failures are due to our seeding of 
the fit which currently tends to favour small-separation geometries. 

Figure 3 shows the region in which the PHOTO star-galaxy separator identified objects 
as galaxies. This figure is of particular interest since, as mentioned above, SDSS quasar 
target selection rejects high-redshift quasar candidates which are classified as galaxies. The 
actual median fiux for high-redshift quasars will be somewhat less than that of the entire 
sample, meaning that fewer resolved objects would be classified as galaxies, but clearly the 
bias against lenses is still quite strong. Figure 4 shows the region in which the PHOTO 
deblends the objects into two separate components. This figure shows the accessible regions 
of the separation - flux ratio parameter space if one wishes to recover wide-separation lenses 
from the SDSS catalogs alone. In these median conditions, even at separations of 3'.'2 the 
PHOTO deblender is only just beginning to deblend point sources with flux ratios less than 
one. Although we concentrated on smaller separation systems. Figure 5 shows, for reference, 
the performance of the PHOTO deblender at larger separations. It is apparent that 4" is 
the approximate splitting at which the PHOTO deblender is complete down to fiux ratios 
which can also be recovered by the Ax^ method. 

Apart from the value of the Ax^ statistic, the best-fit two component model also returns 
a predicted geometry for the system. We can use the known geometry of our simulated 
images to determine the accuracy of the recovered geometry. Figure 6 shows the best-fit 
separation as a function of input separation for all objects with Ax^ > 0.1 (regardless of fiux 
ratio). Apart from a modest systematic tendency to over-estimate the separation of small 
separation pairs, the separation is fit quite accurately. Simply taking the standard deviation 
of the residuals about each point, we estimate the typical error in the separation of be ~ 
0'.'06. Figure 7 shows the logarithmic error in the recovered fiux ratio, again for all objects 
with Ax^ > 0.1 . For each input fiux ratio we fit a Gaussian to the observed distribution of 
errors. Although in each fiux ratio bin there are a few obvious outliers, these Gaussian fits 
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are good enough that we can estimate the errors directly from the measured a. Note that 
the errors in the 1:1 flux ratio bin are asymmetric because the flux ratio is less then one by 
construction. 

We should be careful to note that since the input flux is randomly assigned from a 
Poisson distribution, the purported input flux ratios are only accurate to the level of the 
Poisson noise. Specifically, if T\ and T2 represent the fractions of the total flux in the two 
components, then the uncertainty in the flux ratio from Poisson noise is 



For the median flux, the total counts are ~ 18000, and, for the flux ratios we are 
considering, \JT\T2 ~ 0.5, so that 



Hence the contribution, in quadrature, of the Poisson noise is small enough that we can 
retain our estimate of the error in the logarithm as simply being the values of a shown in 
flgure 7. 

The accuracy of the recovered flux ratios is of vital importance if we wish to use measured 
colour differences between lens components as a selection criterion. We define, for example, 
D{g — r) as the measured g — r colour difference between two lens components. Note that the 
colour difference, D, differs physically from the differential reddening, AE, in that it may 
be caused by intrinsic colour differences, as in the case of a QSO plus star pair, as well by 
reddening. In our notation AD refers to the error in D. Then, if fg and fr are the measured 
flux ratios of two lens components in those bands, we have 




(5) 



or 




(6) 




(7) 




(8) 
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so that the error goes as 



AD{g-r) = 2.5 




(9) 



If we take as a typical value of = cr ~ 0.03, then we get a typical error in the colour 
difference of 



In selecting lens candidates from real data, wc will also want consider the consistency of 
the point source component positions for the best-fit models across the three bands. Even if 
the fiux ratios change between bands, we would still expect the positions of a pair of point 
sources to be the same if they are being successfully recovered in each of the bands. Further, 
it is self-evident that consistent positions are a pre-requisite for the meaningful measurement 
of colours. We define a measure of the difference in component positions between two bands 
as the minimum total distance that the positions of point source components as fit in first 
band must be displaced so as to coincide with the positions of the components as fit in the 
second band, so that 



where the components A and B are chosen so as to minimize the difference. This measure 
considers changes in both separation and position angle, but is independent of changes in 
the flux ratio. In order to measure this difference for our simulated objects, wc produced two 
more sets of simulated pairs of point sources which have identical positions to the set shown 
in figure 2 but which have different reahzations of the Poisson noise. We thus have three 
realizations of each object and we measure the average of the position differences between 
them 



AD{g-r) ~ 0.1 



(10) 



ApositioU]^ 2 = 



(11) 



Aposition = J (position^ 2 + position^ 3 -|- position2 3)78 



(12) 
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Figure 8 shows the distribution of Aposition for those objects in which all three real- 
izations have Ax^ > 0.1 . When selecting objects on consistent geometry, we will require 
Aposition < 0.4 arcseconds. 

Apart from geometry, the total flux of a lens and the seeing conditions under which 
it was observed also greatly affect its detectibility. Figure 9 reproduces the fiducial set of 
simulations shown in figure 2 and illustrates the effect of better and worse seeing as well more 
or less total flux. The seeing values displayed correspond to the 25th-percentile, median, 
and 75th-percentile observed seeing (1'.'8,1'.'6,1'.'4), while the total flux values correspond 
to one magnitude more than the median, the median, and one magnitude less than the 
median. It should be noted that the median observed seeing in the EDR will probably 
end up being good deal worse than the median value for the whole survey. Ultimately, 
for statistical purposes, we envision a suite of simulations which sufficiently-finely samples 
the four-dimensional separation - flux ratio - total flux - seeing parameter space to allow the 
selection function for every quasar in the sample to be accurately determined through simple 
interpolation. 



2.3.2. Four- Component Systems 

While two-component geometries are spanned in the separation - flux ratio parameter 
space, the equivalent parameter space for four-component geometries has eight dimensions, 
albeit with some rotation and reflection symmetries. Consequently, not only is an exhaus- 
tive simulation of four-component geometries computationally-daunting, it would also be 
challenging to interpret the results of such a simulation (eight-dimensional contour plots!). 
Further, as mentioned above, Rusin and Tegmark argue that simple lens models cannot 
reproduce the observed quad fraction for CLASS lenses, meaning that producing a believ- 
able ensemble of lens geometries based on theoretical expectations is at present a dubious 
undertaking. Nonetheless, we wish to at least qualitatively ascertain how well our mor- 
phological selection works for four-component systems. These considerations lead us to the 
following observationally-motivated toy model. We envision a one-dimensional sequence of 
quad-geometries running for highly-symmetric to highly-asymmetric. We chose the known 
lenses H1413+117 (Magain et al. 1988) (symmetric) and PG1115+080 (Weymann et al. 
1980) (asymmetric) as the end-points of our sequence. If we denote a quad as the sum of 
four point-sources in the form 



A I Xi\ 



(13) 



where x and y represent the positions, and F represents the flux ratio of component 
relative to the primary, then we can produce a sequence of quads by interpolating between 
our chosen endpoints so that an interpolated quad will have a geometry given by 

4 / a;,(1413) \ 4 / x,(1115) \ 

Q = E^M ^^(1413) +J](l-a) y,(1115) (0 < a < 1) (14) 

i \ F,(1413) / i \ F,(1115) / 

which is parameterized solely by the interpolation factor a. We used photometry from 
the CASTLES web-site (Kochanek ct al. 2002) to set the flux ratios and relative positions of 
the components of H1413+117 and PG1115+080, but we normalized the total flux of each 
lens to be the same observed median flux we used for our fiducial two-component models. 
Figure 10 illustrates such an interpolated sequence in a pictorial manner. In addition to 
this interpolation, we also allow the scale of the entire system to change, corresponding to 
changing the normalization of the lensing mass. Although the constant size of the pixels 
sets a certain scale, changing the inter-component separations can also be thought of as a 
rough proxy for changing the seeing, since the ratio of the separations to the FWHM changes. 
Figure 1 1 shows the regions of this two dimensional symmetry - size parameter space in which 
our simulated quads pass the Ax^ > 0.1. Here, a scale of one corresponds to a sequence 
of quads where the endpoint geometries have the same size as the real systems H1413+117 
and PG 11 15-1-080. Here we again use the median observed total flux and seeing. Although 
our selection technique is by construction best-suited to detecting two-component objects, 
figure 11 illustrates that it is in fact quite sensitive to any deviations from a point-source 
morphology. Particularly encouraging is the recovery of 1413-type systems which, a priori, 
one might expect to be problematic due to their considerable radial symmetry. 

2.4. Truth Testing 

Formally, the selection function derived above only sets an upper limit on the detectibil- 
ity of a particular lens configuration because our simulations do not include all possibly 
relevant observational effects. For instance, they do not include: i) PSF anisotropy. We use 
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circularly-symmetric PSF's in our simulations whereas the observed SDSS PSFs are often 
(and sometimes highly) non-symmetric, ii) Non-Poisson noise. Bright stars, faint objects, 
and cosmic rays can all contribute sky noise which is not represented by a Poisson distribu- 
tion. Since our minimization is weighted by the square of observed pixel counts, we might 
expect that our fits will be more accurate in simulation where this weighting corresponds 
to the true variance. Even if such effects do not systematically reduce the detectibility of 
non-point source morphologies, they at least are an unmodeled source of noise which might 
statistically reduce the number of lenses detected. 

Another potential difference between our simulations and the real data is the perfor- 
mance of the simplex minimization scheme. We know that sometimes the simplex mini- 
mization will converge to a local, non-global, minimum. For the purposes of measuring the 
selection function, it is irrelevant whether or not this actually happens, only whether or 
not it happens more (or less) frequently for real than for simulated data. In fact, we know 
from the accuracy of the best-fit geometries, figures 6 & 7, that the correct geometry, and 
therefore global minimum, is being found for simulated objects with A^^ > 0.1. Hence, our 
truth testing will tell us whether or not the global minima are being found for the real data 
as well. 

We could test the importance of all remaining unmodeled effects if we were able to 
compare our simulation results to known pairs of point-sources in the SDSS dataset. One 
obvious approach would be to investigate how our two-component fit performs for known 
lenses. Wc performed such a comparison for the newly-discovered lens SDSS1226-0006 (Inada 
et al. 2002a, henceforth 102). This system has two quasar images with a fairly small fiux ratio 
between them. The purpose of this comparison is largely to give the reader a qualitative 
impression of the accuracy of the best-fit geometry . Table 1 shows how the geometries 
recovered from our fit to SDSS imaging compare to those measured from high-resolution 
imaging with the Walter Baade 6.5m telescope as reported in 102. In all three SDSS bands, 
our fit is consistent with the higher-resolution imaging to within the errors derived in section 
2.2 . Table 2 compares the measured Ax^ to a set of simulated objects of the same geometry, 
flux and PSF_W1DTH. The quoted errors for the simulated values are the observed standard 
deviations among 16 realizations of the same object. The simulated values in r and i appear 
to be systematically somewhat higher than the real values, although this may be in part due 
to flux from the lensing galaxy. Unfortunately, few meaningful conclusions can be drawn on 
the basis of a single object. Further, the set of all previously-known lens systems imaged 
by the SDSS to date turns out to be an unsuitable sample for our present purposes; as will 
become evident when we discuss the properties of these objects as observed by the SDSS in 
section 3.3 . 
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As an alternative, we can take advantage of a minor failing of the SDSS galaxy target 
selection to obtain a sample of close pairs of point sources. Since, as demonstrated in the 
preceding section, the PHOTO deblender generally docs not deblend pairs of point sources 
at separations below about 3", binary stars at these separations will be classified by the star- 
galaxy separator as galaxies and, consequently, approximately 1% of all objects targeted in 
the SDSS as galaxies are in fact spectroscopically-identified as stars (Scranton et al. 2002). 
We selected such objects from the EDR and visually-inspected the atlas images to reject 
star/galaxy pairs and other obvious contaminants. This procedure yielded a sample of 307 
objects which, although we cannot be certain of their true geometries, we believe to be at 
least dominated by binary stars. We then fit a two-component model to each of these objects 
to recover a best-fit separation, flux ratio and total flux, as well as a measure of the Ax^. 
We also recorded the PHOTO measured PSF_W1DTH for each object. Next, for each object 
we simulated a pair of point sources with the same scpcration, total flux, and PSF_W1DTH. 
Wc then flt these simulated objects with a two-component model and compared the value 
of the measured Ax^ for the simulated object to that of the real object upon which it was 
modeled. Figure 12 illustrates our truth testing procedure as a flowchart. Figure 13 shows 
the distribution of fractional differences between the simulated and observed value. 



_/'(Ax^) ~ ^^observed ^X-simulated (15) 
^^observed 



If the simulated Ax^ are systematically higher (or lower) than the real values, then 
the center of this distribution should offset from zero. Estimating the center is somewhat 
complicated by the fair number of outhers, particularly for /(Ax^) < -1 (the bins at either 
end of the histogram contain all points at or beyond those values) . We attribute these outhers 
mainly to objects whose underlying geometry is not in reality a pair of point sources. Clearly 
in such a case the Ax^ for the simulated best-flt two-component geometry would be expected 
to be greater than the observed value. However, the goal of this section is to get an estimate 
of the systematic error in Ax^ for cases in which the underlying geometry really is a pair of 
point sources, hence ideally we would like an estimate which ignores such outliers. One way 
to do this is simply to take the median value, which is /(Ax^) = 0.02 . Another is to flt a 
Gaussian to the distribution at —1 < /(Ax^) < 1, as shown in flgure 13. In this case, we 
flnd the center at /(Ax^) = 0.07 . In either case, the systematic error in Ax^ is < 10 %, and 
in fact indicates that the simulations somewhat under-predict the value of Ax^. Ideally, we 
could also compare the dispersion of this distribution to the dispersion of a similiar simulated 
distribution to estimate the random error introduced by unmodelled properties of the data. 
Unfortunately, the objects in this distribution are taken from the main galaxy sample and 
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hence are considerably brighter than the vast majority of the quasars in which we are most 
interested. As a result, for instance, the effect of correctly modelling the sky noise should 
be much less important for these objects. Hence, at this time, we satisfy outselves with 
our conclusion that the systematic differences between the real data and our simulations 
are small (< 10 %) and that, consequently, the selection functions derived in the previous 
sections can be regarded as having been calibrated with real SDSS data. 

3. Application to SDSS Data 

3.1. Observational Details 

The SDSS is a photometric and spectroscopic survey across 10,000 square degrees of the 
northern Galactic cap using the 2.5m SDSS telescope at Apache Point Observatory. SDSS 
imaging is carried out with a wide-field camera (Gunn et al. 1998) with makes nearly simu- 
lataneous observations of objects in five passbands: u g r i z. Together, the passbands cover 
the optical wavelengths from the atmospheric cut-off in the blue to the minimum detectible 
energy for the silicion CCD's in the red (Fukugita et al. 1996). Photometric calibration of 
the imaging survey is seperately carried out by an automated 0.5m telescope (Hogg et al. 

2001) which monitors a set of standard stars (Smith et al. 2002) while photometric data is 
being acquired. The SDSS imaging camera also incorporates astrometric CCDs which pro- 
vide astrometry of detected objects with an acccuracy typically better than O'.'l (Pier et al. 

2002) . 

3.2. Selecting Lens Candidates 

In this section we apply the selection techniques wc described and analyzed in section 
2 to a sample of SDSS quasars. We define as our parent sample the 5120 quasars from the 
EDR imaging which are confirmed by SDSS spectroscopy to be at z > 0.6. We chose to 
restrict ourselves to this sub-sample of the imaging data because at present we are primarily 
interested in presenting methods and the EDR represents a well- vetted test-bed. Our parent 
sample is larger than the 3814 quasars of the SDSS EDR quasar catalog (Schneider et al. 
2002) primarily because we do not restrict ourselves to objects whose spectra appear on 
EDR spectroscopic plates. The SDSS chooses sets of spectroscopic targets (plates) from 
the imaging data (runs). However, not all of the plates targeted from EDR imaging runs 
are included among the plates released as part of the EDR. Hence, our sample additionally 
contains objects whose spectra were targeted on EDR imaging data but not released as 
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part of the EDR. Although this may at first seem a somewhat confusing selection, it arises 
naturally from the fact that we defined our sample solely in terms of the ~ 400 sq. deg. of 
EDR imaging data and then simply found all objects which had SDSS spectroscopy, without 
regard to whether or not those spectra were part of the EDR or not. Another, more minor, 
difference between our sample and the SDSS EDR quasar catalog is that we do not exclude 
any objects based on visual inspections of the spectra or requirements on equivalent widths 
of emission lines. To define a sample of lens candidates, we sequentially apply four cuts to 
our parent sample. Below, we present these cuts together with the number of objects which 
remain after each cut and further comments as appropriate. 

i) Morphological Ax^ Cut: We require that the (S/N)^-weighted average of Ax^ as fit 
in g,r, and i be greater than 0.1. 408 objects remain after this cut. 

ii) Separation Cut: We require that that the average best-fit separation as fit in g,r, 
and i be less than 3'.'0. 105 object remain after this cut. 

As discussed in section 2.1, we concentrated our efforts on recovering lenses with split- 
tings < 3", simply because the vast majority of lenses are expected to have such small 
splittings. Another reason to concentrate on small separation candidates is that while the 
lensing probability falls rapidly with increasing separation, the number of false positives due 
to coincidental alignments rises as the separation squared. Kochanek (1993) calculated in 
detail the expected false positive rate for both stars and galaxies, but here we present the 
following simple estimate: Prom the SDSS number counts of Yasuda et al. (2001), we find 
that at the magnitudes which interest us most (19 < r < 21), the number density of objects 
on the sky is of order 5000 per square deg. Consequently, ~ 1% of quasars will have a 
coincidental companion within 3", a rate which is almost certainly higher than the lensing 
rate. Fortunately, such false positives are less problematic with SDSS data because of the 
colour information available from the five band SDSS photometry (see below), but without 
a seperation cut the efficiency of any selection is certain to be poor. 

The subject of wide-separation candidates requires an understanding of how PHOTO 
deals with multiple objects which are nearby one another on the sky. Here we present a brief 
explanation of parents and children in the SDSS catalogs. For a more complete discussion, the 
reader is referred to Lupton (in preparation). When PHOTO identifies an object, it includes 
an area ^ of the sky around the object in the atlas image. If another object appears within 
that area, then the two (or more) objects are fisted together in the SDSS catalogs as a single 
compound object (parent). Where appropriate, PHOTO then attempts to deblend the parent 
into single objects (children) and, when successful, the deblended children also appear in the 



This area is the predicted full extent of the object based upon the empirical PSF 
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SDSS catalogs. Since a parent object can contain children separated by many arcseconds, 
a large fraction of quasars will have wide-separation parents. SDSS quasar target selection 
only considers the children of deblended parents, but we fit both the children and their 
parent (if any). Consider, for example, a deblended pair of quasars 2" apart; neither child 
is a good lens candidate but obviously the parent would be. Conversely, a non-deblended 
pair of quasars 1" apart may have a star nearby which enters together with the pair into a 
parent object. In this case, both the child (close quasar pair) and the parent (wide quasar(s) 
plus star pair) might pass our morphological selection but only the close pair would pass our 
separation cut. 

iii) Consistent Geometry Cut: We require that Aposition < 0.4 arcseconds geometric 
consistency between the best-fit models in g,r, and i as defined in eqn. 12 . 41 objects 
remain after this cut 

iv) Colour Cut: We require that the colour differences between best-fit model compo- 
nents be consistent with a simple reddening model. Specifically, we require that the measured 
D{g — r) and D{r — i) fall within the region defined by the inequalities 



Note that in order to avoid ambiguities of sign we define colour differences such that 
D{g — r) > and then measure D{r — i) relative to the bluer component in g — r. 13 objects 
remain after this cut. 

We base this cut on the Figure 1 of F99 which indicates that, apart from a pair of 
very highly-reddened systems, the distribution of observed differential extinction for known 
lenses goes to zero at a value of AE(B — T/) ~ 0.4. We devised a simple model to predict 
the colour difference between two quasar images due to reddening in the lens galaxy. In 
what follows, we assume that one of the quasar images is unreddened, so that the colour 
difference is identical to the reddening of the other image. In order to convert between B — V 
and SDSS colours, we need a quasar SED and a reddening law. For redshift < 2, we use 
the composite quasar spectrum of Vanden Berk et al. (2001). Beyond redshift 2, Lyman- 
alpha enters the B-band and this composite does not correctly predict quasar colours due 
to Lyman forest absorption. Instead, at redshift 3, we constructed our own composite using 
the spectra of 198 SDSS quasars with redshifts 2.5 < z < 3.5 . Following F99, we use the 
Galactic reddening law of CardeUi, Clayton, & Mathis (1989) . We chose a value for the 
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i) < D{g-r) + 0.3 
i) > -D{g -r) -0.3 
i) > lD{g-r)- 0.3 
i) < -D{g-r) + 1.3 
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quasar redshift, reddened the composite spectrum as appropriate for a chosen lens galaxy, 
and then convolved the reddened spectrum with the filter transmission curves to find the 
observed colours. By doing this for a number of values of the total extinction, Ay, wc can 
find the values of D{g — r) and D{r — i) which correspond to different values of D{B — V) 
and thus construct reddening vectors in the D{g — r) - D{r — i) colour space. Figure 14 
shows nine reddening vectors of length D{B — V) = 0.5 for quasar redshifts Zqso €{12 3} 
and lens galaxy redshifts Zg G {0.3 0.5 0.7}. The inclusion region described in eqn. 16 is 
also shown. We assume reddening in our galaxy has a constant value of Ay^gaiactic — 0.1 and 
that Rv = 3.1 . F99 found some evidence of differing values of Ry, but presently we do not 
consider this possibility in our simple model. Simply stated, this cut requires that if a colour 
difference is observed in g — r, then the observed colour difference in r — i must be consistent 
with this simple reddening model and that the observed total differential reddening be such 
that D{B -V) < 0.5 

Table 3 summarizes the results of this four-step selection. 

Of the remaining 13 objects, 6 have received follow-up observations which can with 
some certainty confirm or rule out the lensing hypothesis. One of the objects is SDSS1226- 
0006, and five are not lenses. The objects which have been foUowed-up cannot be taken as a 
representative sample because follow-up efforts tend to favour the brighter, wider separation, 
and more equal flux ratio candidates, but if we naively ignore this bias, then wc observe a 
lensing rate of 2 out of 5120, or 4 xlO~^. If the lensing probability follows a Poisson 
distribution, then our sample of 13 candidates implies an upper bound on the observed 
lensing fraction of not more than 22 out of 5120, or 4 xlO~'^ and the one lens in our sample 
imphes a lower bound of not less than 3 xlO~^ (95% confidence levels). 



3.3. Previously-Known Lenses 

An obvious question to ask regarding any selection technique is how well it recovers 
previously-known Icnscd quasars. If very few previously-known systems arc recovered, then 
it stands to reason that some the assumptions we have made in designing our selection tech- 
nique were incorrect. On the other hand, we must guard against overstating the importance 
of recovering known systems. The sample of previously-known systems is very heterogenous; 
subject to a variety of, often unknown, selection effects. Hence, if we were to overly modify 
our selection specifically to recover previously-known systems, we might well introduce those 
same selection effects into our sample. 

At the time of writing, twelve CASTLES lens/binary systems have been imaged by the 
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SDSS. In addition, the current survey area also includes two newly-discovered lens systems, 
SDSS1226-0006, and SDSS0924+0219 (Inada et al. 2002b). Both of these systems were ini- 
tially identified by indcpendant means. Hence, a total of fourteen previously-known systems 
have been imaged at present. Of these, five have image separations greater than 6". Such 
systems are deblended by PHOTO and not in the regime of interest for this work. A further 
three are much too faint (1 magnitude of more below the i = 20.2 flux limit of the high- 
redshift quasar sample) for SDSS spectroscopy and hence out of our parent sample. There 
follows a discussion of the results of our selection for each of the remaining six systems: 

APM08279-|-5255 : (Irwin ct al. 1998) This is a very bright thrcc-imagc lens with 
two image of roughly equal brightness separated by 0'.'38 and a third fainter image (flux 
ratio ~ 11:2) roughly halfway between them (Lewis et al. 2002). This object is actually too 
bright for SDSS spectroscopy and hence is not targeted as a quasar. We recover the correct 
separation and flux ratio of the two brighter components in r and i, but not in g where the 
object is almost two magnitudes fainter {g,r,i — 17.2, 15.4, 14.9) and the seeing is worse. 
Consequently, it passes our Ax^ and consistent geometry cuts but fails the colour cut due to 
a large, best-flt, g — r colour difference. Given the image seperation, our ability to recover 
the geometry at all is quite encouraging, albeit at extremely high signal-to-noise. 

SBS0909-|-532 : (Kochanck ct al. 1997) This is a bright two-image lens with an image 
separation of I'.'l. It is targeted by SDSS quasar target selection and passes all of our lens 
candidate selection criteria. 

SDSS0924-I-0219 : This unusual system appears to have only three bright quasar 
images, but is in fact almost certainly a four image system with the fourth image obscured 
and/or demagnifled by micro-lensing. It is targeted by SDSS quasar target selection and 
passes all of our lens candidate selection criteria. This object demonstrates that our method 
can recover lenses with more complicated geometries than a simple two-image conflguration. 

SDSS1226-0006 : This is a two-image lens with an image separation of 1'.'2. It is 
targeted by SDSS quasar target selection and passes all of our lens candidate selection 
criteria. 

SBS1520-I-530 : (Chavushyan et al. 1997) This system has a rather complex morphol- 
ogy (there are up to four quasar images and up to three lens galaxies). Furthermore, it is 
within 14" of a bright (r = 13.6) star. Presumably due to being in the PSF wing of this star, 
the object is badly- deblended into four children of non-astrophysical colours. Consequently, 
the results of our fltting to this object are not meaningful. The object is also not targeted 
by SDSS quasar target selection. 

HST15433-|-5352 : (Ratnatunga et al. 1999) This is a faint two-image lens system 
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with an image separation of l'/2. The lensing nature of this object is still not certain. It is 
slightly too faint for SDSS quasar target selection {i — 20.6). Additionally, it is classified 
by the PHOTO star-galaxy separator as a galaxy, and hence, it would not be targeted as 
a high-redshift quasar due to its extended morphology. This is an example of the anti-lens 
bias of the high-redshift quasar sample which was discussed in section 2.3.1 . This object 
passes our Ax^ and consistent geometry cuts but fails the colour cut due to a large best-fit 
r — i colour difference. This colour difference is due to a somewhat inaccurate fit in r which 
is caused by the lower signal-to- noise in this band (r = 21.6). 

The results of our selection for these six objects are summarized in table 4. For the 
three objects which we do select, our selection techniques work just as envisioned; we obtain 
a robust Ax^ in all three fitted bands and the measured colour differences are small. For the 
two objects which are not selected, inaccuracies in the best-fit model in one or more bands 
lead to large colour differences which cause the object to be rejected as a lens candidate. 
These two examples show that our colour cuts are in fact also cuts on the accuracy of the 
two-component model when the object lies near the selection threshold. However, given that 
one of these two failures occurs at a very small splitting angle and the other at a signal-to- 
noise below the limit of the SDSS quasar sample, we feel that neither warrants changing our 
selection criteria. That said, as more known systems are imaged, we intend to repeat this 
exercise and re-appraise our selection techniques accordingly. 



4. Discussion 

The main goal of this work has been to present our method for quantifying the selec- 
tion function of our (or any) selection algorithm. We have explicitly presented the selection 
function only for a somewhat coarse sampling of the available parameter space, particularly 
for four-component systems, but our simulation method can be straight-fowardly applied to 
determine the selectibihty of any lens geometry. The accuracy to which the selection can 
be determined is limited by the behaviour of unsimulated properties on the real data, as 
considered in section 2.4 . It ceases to be useful to further improve one's understanding of 
the selection function once the error associated with it is much less than the uncertainty 
introduced in the statistics of the lens sample by Poisson fluctuations. If the lensing rate is 
indeed < 10^'\ then the Poisson uncertainty for the complete sample of 10^ SDSS spectro- 
scopic quasars will be > 10 %. Such an accurary should be readily-accessible through our 
methods, at least for the morphological selection. 

In order to make use of our selection function, one needs to predict the distribution 
of lens geometries so as to correct the observed distribution for those lenses which will not 
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be selected. We have not undertaken such an exercise at present, consequently our upper 
bound on the observed lensing fraction is only a lower bound on the upper bound on the true 
lensing fraction. Presently, we feel that attempting to compute the selection function for the 
whole survey would be premature. There are still a few obvious improvements to our selection 
procedure, such as an improved seeding scheme to recover wider separation candidates, which 
should be made. Also, as stated above, the data, in particular the seeing conditions, of the 
EDR are probably not representative of the whole survey. Additionally, since the prediction 
of lensing statistics requires knowledge of the quasar luminosity function, it is sensible to 
await determination of the SDSS luminosity function since, almost by definition, it will be 
luminosity function best-suited to the this aim. 

A theoretical model of the lens galaxy population will also be required to better constrain 
the effects of the lens galaxy flux upon the selection function. The estimates in section 
2.1 imply that our simplifying assumption that the lens galaxy flux is negligible probably 
introduces moderate errors in the selection function for the majority of quasars at redshifts 
less than 3, and significant ones at higher redshifts. A more detailed calculation of this 
effect is required. If, in fact, many lenses have prominent lens galaxies, then our simulation 
techniques could still be used to measure the selection function, but the accuracy of the 
result would be dependant on the accuracy with which we could predict the photometric 
properties of the lens galaxies. 

Another unresolved problem is the use of colour as a selection criterion. On the one 
hand, to not use the multiple colours available from SDSS imaging is a enormous waste of 
information, particularly when attempting to reject stellar interlopers. For that matter, we 
have been conservative in this work by not using u-hand imaging although, for quasars with 
z < 2.5, u — g colour is the strongest discriminant between quasars and stars. On the other 
hand, while our simulations allow for an accurate quantification of the selection function 
associated with our morphological criterion, the selection effects associated with colour se- 
lection remain difficult to similarly characterize. Further, if, contrary to our expectations, a 
sizable fraction of lenses are in fact being reddened, or otherwise displaced in colour space, 
so as to fall onto the stellar locus, and hence out of the quasar sample, then there is little 
prospect of ever accurately measuring the lensing rate. 

There are a number of obvious observational extensions to the work presented in this 
paper. First would be to complete the follow-up of the candidates identified in section 3 so as 
to get an actual measure of the observed lensing rate for this sample. Beyond that, applying 
our selection techniques to the full set of available SDSS quasars (there are approximately 
30000 at the time of writing) would produce a sample of candidates which should give a 
more authoritative estimate of the lensing rate. A sample of this size, hkely several hun- 
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dred candidates, would probably best be approached by beginning with a program of high- 
resolution imaging to confirm candidate geometries before attempting the telescope-intensive 
high spatial-resolution spectroscopy or near-IR imaging which is required to definitively test 
the lensing hypothesis. 

As shown in section 2.3, many of the most prominent lenses at redshifts greater than 3 
will not be in the SDSS quasars sample due to their extended morphology. However, due to 
higher lensing rate at these redshifts, it may be profitable to attempt to identify such lenses 
from the SDSS imaging alone. In particular, if our morphological criterion can be refined to 
allow for the selection of lens geometries from among normal, red galaxy morphologies, then 
it may be possible to carry out an cfiicicnt survey for lenses at these redshifts. 

We draw attention to the possibility of the application of our simulation techniques to 
problems beyond the discovery of gravitational lenses. Other researchers have previously 
used simulated SDSS images to study selection effects (cf. Kim et al. 2002), but the sim- 
ulation procedure we have described has the advantage of allowing one to create simulated 
images of any sort of object with realistic noise properties, particularly when the simulations 
are calibrated by real data as in section 2.4 . Such calibrated images are well-suited for 
any undertaking which involves additional processing of the SDSS imaging data. Bulge-disk 
decomposition and other measures of galaxy shapes beyond the measured PHOTO parame- 
ters are obvious examples of problems which might be approached in this way. Additionally, 
because the simulations are passed through the same reductions as the real data, one has 
the ability to investigate how a given class of objects is treated by PHOTO. One possibly 
interesting application of this would be to investigate the treatment of bright galaxies by the 
PHOTO deblender. 

Comparison of our preliminary observed lensing fraction with previous results is com- 
plicated by the absence of a comparable homogeneous sample. The HST Lens Snapshot 
Survey (Maoz et al. 1993) found 6 lenses among 502 optically-selected quasars (a rate of 1.2 
xl0~^), but the magnification bias is likely much more pronounced for their brighter parent 
sample. The CLASS survey found 22 lenses among 11,685 fiat-spectrum radio sources (a 
rate of 1.9x10"^), but in this case the uncertain redshift distribution and optical luminosity 
function of their soTirccs again makes direct comparison difficult. A definitive measurement 
of the lensing rate from the SDSS will not only await further survey data, but also a great 
deal of observational follow-up outside of the scope of the survey itself. 
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A. Simulating SDSS Images 

In this appendix, we describe how we created the simulated SDSS images which we used 
in section 2 to measure the selection function of our morphological selection criterion. We 
begin by choosing a value for the PSF FWHM and creating a simulated image of a field of 
bright point sources with that FWHM. The image is generated by first integrating analytic 
PSFs on a grid corresponding to the SDSS pixel scale of Of! A /pixel. We used an analytic 
PSF composed of two Moffat functions with Pi — 7 and = 2, as recommend by Racine 
(1996). It should be noted that PHOTO uses a double Gaussian model when estimating 
the PSF_WIDTH in a frame and, consequently, we find that the measured PSF_WIDTH is 
systematically greater than the input FWHM of our Moffat functions. To first order, the 
relation between the two is PSF_WIDTH ~ 1.15 x FHWMMoffat- Throughout the previous 
sections of this paper quoted FWHM values have been for the PSF_WIDTH as reported by 
PHOTO. We then add a constant sky level and assign 'observed' counts for each pixel as a 
random Poisson variable with a mean value of the object plus sky counts. For the PSF fields 
we used the observed median EDR sky value in r. We also add a Poisson random variable 
corresponding to the observed dark variance. Having recorded our field of PSF stars as a 
standard FITS file, we then pass this simulated field through to the so-called fake stamp 
collection pipeline (fsc_pipeline) which produces ancillary files required for running PHOTO. 
The image is then passed through the postage stamp pipeline (ps_pipeline) which outputs an 
empirical PSF (psField). It should be noted that, in contrast to the PSF_WIDTH estimate, 
this empirical PSF does not assume a double Gaussian or any other particular functional 
form (again, see Lupton et al. 2001 for details). The procedure for simulated lens candidates 
is the same, except that instead of single point sources we lay down a pair of point sources 
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with the desired separation, flux ratio, and signal-to-noise. Also, for simulated ensembles 
where each object is reproduced in ten iterations, we folded-in the observed distribution 
of sky values by assigning to the first iteration, the 5th-percentile sky level; to the 2nd 
iteration, the 15th-pcrcentile sky level, and so on. The sky distribution we used was simply 
the observed sky values in the 5120 fields which contain our SDSS quasar sample. Again, 
we prepare to run PHOTO (fsc_pipeline) and then we pass our image through the frames 
pipeline (frames_pipeline) which, making use of the empirical PSF which we just previously 
constructed, outputs the cropped atlas images (fpAtlas) as well as the measured PHOTO 
outputs, including the results of the deblending routine. In regular survey operations, the 
objects and PSF field stars are in the same frame, but we perform these steps separately for 
efficiency as the PSF need only be measured once per choice of the FWHM. Note that we 
do not complete the standard photometric pipeline by calibrating the objects (nfcalib) as 
photometric calibration is not of interest to us. Figure 15 illustrates the simulation procedure 
as a fiowchart. The fiowchart also includes our psf- fitting to the atlas images. 
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Fig. 1. — Solid histogram shows the distribution A^;^ for 5120 simulated point sources with 
flux, seeing, and noise chosen to match the sample of SDSS EDR quasars. The dashed line 
is a Gaussian-fit to this distribution used to estimate the false positive rate as a function of 
the Ax^ threshold. 
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Fig. 2. — Regions of the separation - flux ratio parameter space in which the morphological 
Ax^ > 0.1 cut recovers pairs of point sources for the median seeing and flux of SDSS EDR 
quasars. Solid red circles indicate points where out of 10 iterations are recovered. Open 
green stars indicate points in which between 1 and 9 out of 10 iterations are recovered. Open 
blue circles indicate points in which 10 out of 10 iterations are recovered. 
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Fig. 3. — Same as Figure 2 except that the PHOTO star-galaxy separation is used as a 
selection criterion instead of the Ax^ statistic. Selected points (open circles and stars) are 
geometries classified as galaxies by PHOTO 
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Fig. 4. — Same as Figure 2 except that the PHOTO deblender is used as a selection criterion 
instead of the Ax^ statistic. Selected points (open circles and stars) are geometries deblended 
by PHOTO 
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Fig. 5. — Same as Figure 4 except that the x-axis (separation) is expanded by a factor of 
two. 
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Fig. 6. — Upper panel shows the recovered separation of simulated pairs of point sources 
as a function of the input separation. Solid circles indicate the input separation. The open 
circles indicate the average recovered separation of all objects at a given input separation 
which have Ax^ > 0.1 . The error bars are the measured standard deviation. Lower panel 
shows the difference between the input and recovered separations. 
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Fig. 7. — Distribution of logarithmic errors in the recovered flux ratio for all simulated pairs 
of point sources which have Ax^ > 0.1 . The four panels show the distributions for different 
values of the input flux ratio as labeled on the y-axis. The dashed curve shows a Gaussian 
flt to each of these distributions. The a for each Gaussian is also indicated. 
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Fig. 8. — Distribution of Aposition for three- fold realizations of pairs of simulated point 
sources where all three reahzations have Ax^ > 0.1 . 
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Fig. 9. — Effect of changing the seeing and total flux upon the selection function. The center 
panel is identical to Figure 2. Going from left to right, the flux decreases by one magnitude 
(left), is the median (center), increases by one magnitude (right). Going from the bottom to 
top, the seeing is the 25-th percentile value (I'.'S) (bottom), is the median (I'.'G) (center), is 
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Fig. 10. — Cartoon illustrating a series of five interpolated four-component lens geometries 
running from less symmetric (PG1115-|-080,a=0) to more symmetric (H1413-|-117,a=l). The 
radii of the circles indicates the relative fluxes of the components. 
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Fig. 11. — Regions of the symmetry - scale parameter space in which the morphological 
Ax^ > 0.1 cut recovers interpolated quads for the median seeing and flux of SDSS EDR 
quasars. Along the x-axis the scale of the system is changed between 0.2 and 2.0 of the actual 
observed size. Along the y-axis the interpolation factor is changed from (asymmetric 1115) 
to 1 (symmetric 1413). The meaning of the points is the same as in Figure 2. 



(K! 



to 



o 

I 

O 

i-i 

O 

a> 
o 



2 

o 



^ 8 

^ d 



SDSS Image of 

Point-Source 

Pair 




5' 

CI- >Td 

p i-i 

rt- O 



I-i 

o 
o 



I-i 



Oi 



Measured 
Geometry and 
Delta Chi-Squared 


—> 


Comparison of 
Delta Chi-Squared 




/ 




/ 


\ 


Simulated Image 
with same 
Geometry, Seeing, 
Signal-to-Noise 






/ 






Two-Component 
Fit 


—> 


Measured 
Geometry and 
Delta Chi-Squared 



"T 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 r 



30 



w 
o 

o 



20 



o 

Dh 



CO 

S-H 
•rH 

CO 
O 

Sh 

0) 
X! 



_l I I I I I I I I I I I L 



-0.5 







J I I L_ 



0.5 



Fig. 13. — Histogram of the distribution of fractional differences in tlie measured value of 
Ax^ between simulated images and real SDSS imaging data. The dashed line is a Gaussian 
fit to the distribution for values — 1 < /(Ax^) < 1. 
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Fig. 14. — Reddening vectors of length AE{B — V) — 0.5 for the composite quasar spectrum 
of Vanden Berk et al. (2001) with quasar redshifts Zgso €{12 3} and lens galaxy redshifts 
Zg e {0.3 0.5 0.7} (solid curves). The dashed lines indicate the inclusion region demarked 
by colour cuts described in eqn. 16 . The circle in the lower right has radius of 0.1, equal to 
the estimated error in the colour difference measurements, as derived in eqn. 10 . 
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Fig. 15. — Flowchart illustrating the steps in our simulation and psf-fitting procedures. The 
simulation inputs are in the uppermost box and the psf-fitting outputs in the lower most 
box. 
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Table 1. Comparison of Model to Data for SDSS1226-0006 





Band 


Separation 
(") 


Flux Ratio 


Observed'^ 


9 


1.24 ± 0.02 


1.8 




i 


1.24 ± 0.02 


2.2 




9 


1.28 


1.7 


Fit to SDSS Data 


r 


1.24 


1.8 




i 


1.21 


2.1 



'^Walter Baade 6.5m with MagIC 
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Table 2. Comparison of Simulation to Data for SDSS1226-0006 





Band 


Ax' 




9 


7.42 


Fit to Real SDSS data 


r 


8.34 




i 


5.31 




9 


7.64 ± 0.30 


Fit to Simulated SDSS data 


r 


10.55 ± 0.35 




i 


7.01 ± 0.33 
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Table 3. Summary of Selection Procedure 



Select On 


Number of 
Remaining Candidates 


Percentage of 
Candidates Remaining 


Spectrum (QSO z > 0.6) 


5120 


100 % 


Ax^ > 0.1 


408 


8% 


Separation < 3" 


103 


2% 


Consistent Geometry Aposition < 0.4 


41 


0.8 % 


Consistent Colour 


13 


0.3 % 
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Table 4. Selection Results for Previously-Known Lenses 



Name 


Selected by SDSS 
quasar target selection 


Selected as 
Lens Candidate 


APM08279+5255 


N 


N 


SBS0909+532 


Y 


Y 


SDSS0924+0219 


Y 


Y 


SDSS1226-0006 


Y 


Y 


SBS1520+530 


N 


N 


HST15433+5352 


N 


N 



